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ABSTRACT 

A  procedure  has  been  developed  for  treating  the  dynamic 
interaction  between  a  separated  wake  and  a  surface  undergoing  an 
unsteady  motion.  The  basis  of  the  method  is  an  unsteady  (time¬ 
stepping)  panel  method  coupled  with  unsteady  integral  boundary 
layer  codes.  Pilot  codes  have  been  developed  for  both  two-  and 
three-dimensional  conditions.  Results  presented  here  are  mainly 
from  the  two-dimensional  code  in  which  the  various  routines  for 
controlling  the  dynamic  wake  model  have  been  developed.  Some 
viscous/inviscid  three-dimensional  results  are  shown.  The  long 
term  objective  is  to  treat  complete  aircraft  configurations 
through  high  angle-of-attack  maneuvers. 


1.0  INTRODUCTION 


Recent  interest  in  the  concept  of  supermaneuverability 
( e.g (1)  and  (2))  has  stimulated  renewed  efforts  to  explore 
ways  of  predicting  the  unsteady  flow  characteristics  of  combat 
aircraft  in  high  angle-of-attack  maneuvers.  The  problem  requires 
treatment  of  general  configurations  having  separated  flows  whose 
extent  varies  with  time  and  whose  downstream  development  cul¬ 
minates  in  complex,  energetic  vortical  structures  which  strongly 
interact  with  the  aircraft  surfaces.  It  is  this  interaction 
which  greatly  affects  the  dynamic  behavior  of  the  aircraft  and 
which  needs  to  be  understood  and  to  be  predictable  if  the  super- 
maneuverable  concept  is  to  be  exploited  and  integrated  into  new 
designs. 

Strictly  the  above  problem  requires  treatment  by  the  Navier- 
Stokes  equations,  but  such  solutions  are  very  expensive  and  are 
currently  limited  to  fairly  simple  configurations  (e.g.,  see  the 
recent  review  by  Laschka  (3)).  The  geometric  complexity  of 
modern  combat  aircraft  is  clearly  beyond  practical  treatment  by 
Navier-Stokes  solvers  and  is  likely  to  remain  so  for  some  time  to 
come.  The  present  project,  therefore,  examines  an  alternative 
approach  based  on  coupling  between  unsteady  boundary  layer 
methods  and  an  unsteady  (time-stepping)  panel  method  which  in¬ 
cludes  a  dynamic  wake  model. 

Panel  codes  such  as  VSAERO  (4)  routinely  handle  the  geometry 
of  complex  aircraft  configurations.  The  calculated  steady  and 
quasi-steady  aerodynamic  characteristics  include  corrections  for 
■real  flow"  effects:  iteratively  coupled  boundary  layer  codes 
and  a  wake  relaxation  routine  modify  the  basic  inviscid  solution 
and  the  effects  of  compressibility  are  treated  by  approximate 
corrections  such  as  those  of  Prandtl-Glauert  or  Karman-Tsien. 
Also,  separated  flow  problems  have  been  successfully  treated 
using  vortex  sheet  modeling  of  the  free  shear  layers  "enclosing" 
the  separated  regions  (5).  The  low  computing  cost  of  the  basic 
VSAERO  method  has  allowed  practical  time-stepping  calculations  to 
be  performed  with  fairly  complex  geometries  (6). 

The  separated  flow  and  unsteady  capabilities  have  now  been 
combined  in  a  new  formulation  installed  in  both  two-  and  three- 
dimensional  pilot  codes.  These  have  already  demonstrated  a  cost 
effective,  powerful  capability  (7).  The  present  approach  avoids 
the  empiricisms  of  many  of  the  earlier  methods  for  unsteady 
separated  flows  (e.g.,  see  the  reviews  by  McCroskey  (8),  (9)). 
It  extends  the  unsteady  capability  for  thin  wings  offered  by 
various  vortex-lattice  methods,  e.g.,  Rehbach  (10),  Kandil  (11), 
and  Levin  and  Katz  (12),  to  the  case  of  thick  wings  with  fuse¬ 
lages,  etc.,  comprising  a  complete  aircraft  configuration.  The 
greater  details  of  surface  pressure  and  velocity  distributions 
offered  by  the  approach  allow  viscous  effects  to  be  included 
through  coupled  unsteady  boundary  layer  methods.  In  this  respect 
the  present  approach  parallels  the  two-dimensional  work  by 
Strickland  et  al.  (13)  and  the  investigations  by  Geissler  (14). 


The  long  term  objective  of  the  present  approach  is  to  pro¬ 
vide  a  practical  engineering  tool  capable  of  treating  complete 
aircraft  configurations  in  arbitrary  motion.  While  separation 
lines  are  provided  by  the  unsteady  boundary  layer  calculations, 
the  approach  is  concerned  more  with  the  effects  of  the  large- 
scale  vortical  structures  and  their  dynamic  interactions  with  the 
vehicle  surfaces  than  with  the  small-scale  details  of  the  viscous 
regions.  However,  the  method  should  provide  a  basis  for  future 
couplings  with  more  complete  treatment  of  the  surface  viscous 
effects  of  unsteady  separating  flows  such  as  those  being  pursued 
by  Geissler  (14),  Cebeci  et  al.  (15)  and  Le  Balleur  (16),  among 
others.  In  addition,  the  method  provides  a  foundation  for  ap¬ 
plying  improved  modeling  techniques  resulting  from  further 
evaluation  of  experimental  results  such  as  those  of  McCroskey  and 
Philippe  (17),  McAlister  and  Carr  (18),  Francis,  Keese  and 
Retelle  (19),  Helin  and  Walker  (20),  and  Graham  and  Strickland 
(21) . 

Finally,  the  approach  presented  here  has  a  potential  appli¬ 
cation  to  a  broad  range  of  unsteady  problems  involving  arbitrary 
geometries  and  arbitrary  motions. 


2.0  METHOD  FORMULATION 


The  problem  is  first  reduced  to  a  flow  model  in  which  the 
aircraft  and  its  wake  are  immersed  in  an  inviscid,  irrotational 
and  incompressible  flow;  regions  dominated  by  viscous  effects  are 
assumed  to  be  confined  to  the  wake  and  to  thin  boundary  layers  on 
the  aircraft  surfaces.  Furthermore,  the  separated  wakes  are  as¬ 
sumed  to  be  contained  within  thin  shear  layers — represented  by 
vortex  sheets.  These  sheets  are  allowed  to  distort  and  to  roll 
up  to  form  discrete  vortices  which,  for  numerical  purposes,  are 
given  a  small  but  finite  "core".  The  wake  fluid  between  the 
vortex  sheets  and  external  to  the  vortex  cores  is  assumed  to  be  a 
low-energy  potential  flow.  Figure  1  illustrates  the  features  of 
the  model.  It  shows  a  section  taken  through  a  wing  undergoing  a 
high  a  pitch  maneuver.  The  wake  vortex  sheets  contain  the  start 
of  vortex  roll-up  due  to  the  transient  in  the  rate  of  vorticity 
shedding.  This  simple  model  of  an  extensive  separated  region  has 
been  used  successfully  in  the  "steady"  separated  flow  case  (5), 
(22),  and  for  unsteady  separated  flows  from  thin  sharp-edge 
surfaces  (12).  It  is  justifiable  based  on  experimental  observa¬ 
tions  that  the  rate  of  vorticity  dissipation  in  the  separated 
free-shear  layers  is  initially  very  low. 

Consider  the  external  flow;  the  initial  assumptions  above 
allow  a  velocity  potential,  <J> ,  to  be  defined,  and  the  continuity 
equation  reduces  to  Laplace's  equation. 

V2$  =  0 

This,  coupled  with  Green's  Theorem,  yields  the  basic  equa¬ 
tion  for  the  panel  method  (e.g.,  see  (6)).  The  form  used  here 
employs  the  internal  Dirichlet  boundary  condition  of  zero  pertur¬ 
bation  potential  applied  at  points  P  on  the  inside  of  the  air¬ 
craft  surface,  S: 
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V(^)dW  -  0 
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where  S-P  signifies  that  the  integral  excludes  the  point,  P, 
whose  contribution  is  given  by  the  second  term,  n  is  the  unit 
normal  to  the  surface,  S,  or  to  the  wake,  W;  r  is  the  distance 
form  an  inducing  element,  dS  or  dW,  to  the  point  P;  u  is  the 
surface  perturbation  potential  divided  by  4t;  and  uy  is  the  jump 
in  tiLLfll  potential  across  the  wake  sheet  divided  by  4~.  o  ,  the 
source  strength,  is  equal  to  the  normal  component  of  the  pertur¬ 
bation  velocity,  v,  divided  by  4^;  i.e., 
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Figure  1.  Features  of  the  Simplified  Flow  Model. 


It  is  observed  from  Eq.  (1)  that  in  the  present  formulation 
the  wake  is  represented  by  doublet  singularities  only  (the  gra¬ 
dient  of  which  gives  the  velocity  jump  across  the  wake  sheet  and, 
hence,  the  local  vorticity  vector)  and  the  aircraft  surface  is 
represented  by  doublet  and  source  singularities.  The  basic  un¬ 
known  is  the  surface  doublet  distribution.  The  source  distribu¬ 
tion  can  be  evaluated  directly  from  the  external  Neuman  boundary 
condition  applied  at  each  instant  of  time  as  described  below. 

The  aircraft-fixed  reference  frame  (x,y,z)  is  moving  and 
rotating  relative  to  an  inertial  frame  (X,Y,Z),  Figure  2. 

The  flow  velocity  relative  to  the  aircraft-fixed  frame  is 

V  .  v  -  VA  -  &h  -  r 

where  VA  is  the  aircraft  frame  velocity  viewed  from  the  inertial 
frame.  ft  is  the  angular  velocity  about  a  unit  axis,  h,  passing 
through  the  aircraft  frame  origin,  and  r  is  the  position  vector 
of  the  point  in  the  aircraft  reference  frame  where  the  velocity 
is  being  evaluated. 

The  normal  component  of  V  at  a  point  on  the  surface  is, 
therefore: 

VN  «  n  .  V 

*  n  •  v  -  n  •  VA  -  ftn  •  h~r 
and  using  Eq.  (2),  the  source  value  becomes 

0  *  ^T|vN  ♦  n  •  VA  ♦  Qn  •  h  ~  rj  (4) 

The  resultant  normal  flow  velocity,  Vn,  relative  to  the 
aircraft  surface  is  generally  zero.  However,  Vn  can  take  non¬ 
zero  values  to  represent  regions  of  outflow  (Vn  >  0)  or  regions 
of  inflow  (Vn  <  0)  to  model  exhaust  and  inlet  flows,  respective¬ 
ly;  or  to  represent  boundary  layer  displacement  effects;  i.e., 

V  =  3/9s  (V  6*) 

N  e 

where  Ve ,  6*  are  the  local  external  velocity  magnitude  and 

boundary  layer  displacement  thickness,  respectively.  The  deriva¬ 
tive  is  taken  with  respect  to  distance,  s,  along  the  external 
streamline.  The  surface  distribution  of  V^  is  provided  by  the 
boundary  layer  calculation  in  the  previous  time  step. 

When  a  doublet  solution  has  been  obtained,  velocities  are 
calculated  at  points  on  the  wake  vortex  sheets  using  Eq.  (3),  and 
the  complete  wake  is  convected  downstream  through  a  small  time 
step.  A  new  strip  of  wake  is  formed  along  the  separation  lines 
as  described  in  the  next  section. 


INERTIAL  FRAME  AIRCRAFT  REFERENCE  FRAME 


Figure  2.  Aircraft  Movement  Relative  to  the  Inertial  Frame. 
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3.0  NUMERICAL  PROCEDURE 

The  method  is  divided  into  a  number  of  steps  as  shown  in 
Figure  3.  The  initial  setup,  geometry  definition  and  formation 
of  the  basic  matrix  of  influence  coefficients  for  the  aircraft 
surface  occur  outside  the  time  step  loop.  The  calculations  can 
start  either  with  a  wake,  i.e.,  an  initial  steady-state  condi¬ 
tion,  or  without  a  wake.  The  latter  implies  an  impulsive  start. 
The  calculations  then  proceed  through  the  potential  flow  solution 
and  boundary  layer  analysis.  The  new  wake  is  then  formed  while 
taking  a  small  time  step.  Wake  influence  coefficients  are  then 
j  calculated  and  combined  with  the  basic  matrix  to  generate  the 

1  next  potential  flow  solution.  The  subsequent  discussion  outlines 

*  these  steps  in  more  detail. 

3.1  Initial  Setup 

The  initial  setup  requires  user  input  to  establish  the 
general  conditions  for  the  calculations,  e.g.,  reference  quanti¬ 
ties,  type  of  problems,  form  of  output  required,  etc.  In  par¬ 
ticular,  the  motion  and  orientation  of  the  aircraft  frame  rela¬ 
tive  to  the  inertial  frame  is  established  here  as  a  function  of 
time.  This  is  assumed  to  be  in  the  form  of  a  schedule  tabulated 
with  time.  For  simple  motions  or  for  completely  arbitrary  mo¬ 
tions,  the  schedule  can  be  input  directly,  e.g.,  a  schedule  of 
time  and  a  provides  a  pitch  motion.  Tables  for  certain  special 
motions,  e.g.,  sinusoidal  oscillations  coning,  etc.  may  be 
generated  by  a  subroutine  given  certain  key  parameters.  Inside 
the  program  the  time  derivatives  are  evaluated  using  central 
differencing  in  the  table. 

3.2  Geometry  Definition 

The  geometry  of  the  aircraft  is  specified  next.  This  is 
described  within  a  geometry  reference  frame.  Figure  4,  which  is 
selected  by  the  user  for  convenience.  This  is  a  body-fixed  frame 
and  may  be  parallel  to  but  shifted  from  the  aircraft  frame  which 
normally  has  its  origin  at  the  center  of  gravity.  The  user 
supplies  the  transformation  (usually  xcg^  to  move  from  the  geo¬ 
metry  frame  to  the  aircraft  frame. 


Once  the  surface  geometry  is  defined,  the  program  generates 
the  quadrilateral  panels  which  represent  the  aircraft  surface. 


3.3  Matrix  of  Influence  Coefficients 

Each  surface  panel  is  assumed  to  have  a  uniform  source  and 
doublet  distribution.  The  integrals  in  Eq.  (1)  can  then  be 
evaluated  in  closed  form  for  each  panel.  A  control  point  where 
Eq.  (1)  is  applied  is  located  under  the  center  of  each  panel.  If 
there  are  N  panels  representing  the  aircraft  surface,  Eq.  (1) 
becomes 
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Nw  is  the  number  of  wake  panels  at  this  step.  Cjk ,  Bjk  are  the 
influence  coefficients  for  unit  doublet  and  source  uniform  dis¬ 
tributions  on  a  quadrilateral  panel.  oK  is  the  instantaneous 
panel  source  value  from  Eq.  (4).  This  is  regarded  in  several 
parts  for  convenience  in  the  coding. 

Wake  panels  have  a  constant  doublet  distribution  in  the 
normal  to  streamwise  direction  and  a  linear  doublet  distribution 
(i.e.,  constant  vorticity)  in  the  streamwise  direction.  Hence, 
the  wake  integral  term  in  the  quantity,  Ej,  in  Eq.  (5)  has  an 
additional  part:  ywk  is  the  vorticity  value  on  the  Kth  wake 
panel  and  DjK  is  the  vorticity  influence  coefficient  for  the 
quadrilateral  panel.  The  vorticity  vector  is  normal  to  the  local 
stream  direction;  its  modulus  is  the  local  negative  gradient  of 
the  doublet  distribution  down  the  center  of  each  column  of  wake 
panels. 

The  doublet  values  on  the  wake  panels  are  known  from  solu¬ 
tions  at  earlier  time  steps.  Only  the  values  along  the  separa¬ 
tion  lines  (including  trailing  edges)  are  unknown  at  each  time 
step  and  these  are  obtained  from  the  surface  doublet  solution. 
Actually,  the  values  are  taken  from  the  local  jump  in  total 
potential  between  the  pair  of  wake  shedding  panels  to  either  side 
of  the  separation  line. 

3.4  Potential  Flow  Solution 

The  solution  to  Eq.  (5)  provides  the  surface  doublet  values. 
The  surface  velocities  and  pressures  can  then  be  evaluated.  The 
incompressible  pressure  coefficient  is 


CP  -  (v*  -  Va  ♦  2  /*fc)/V* 
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where  Vg  *  +  fih  ~  r  is  the  instantaneous  velocity  of  a  point 
on  the  surface  relative  to  the  inertial  frame  and  V  is  the  modu¬ 
lus  of  the  flow  velocity  from  Eq.  (3).  The  perturbation  veloci¬ 
ty,  v,  is  evaluated  from  the  gradient  of  u  at  each  panel.  A 
local  second-order  fit  to  the  panel  doublet  values  in  two  direc¬ 
tions  is  differentiated  to  obtain  this. 

The  pressure  distribution  can  be  integrated  over  the  air¬ 
craft  surface  to  yield  the  instantaneous  force  and  moment  coeffi¬ 
cients.  These  are  saved  in  a  file  for  the  force  and  moment 
history. 

Next,  the  surface  velocity  distribution  is  examined  to  lo¬ 
cate  lines  of  attachment  (divergence)  and  lines  of  convergence 
(separation)  in  the  external  flow.  In  the  three-dimensional 
version,  surface  streamlines  (external  to  the  boundary  layer)  are 
also  computed.  With  this  information  sets  of  points  and  their 
velocities  are  assembled  for  the  boundary  layer  analysis. 


Curie's  original  method  (23)  has  been  modified  to  calculate 
the  unsteady  boundary  layer  development.  This  is  achieved  by 
solving  the  unsteady  momentum  integral  equation  using  a  Runge- 
Kutta  method.  The  turbulent  boundary  layer  method  is  based  on 
the  unsteady  momentum  integral  equation  as  in  the  laminar  bound¬ 
ary  layer  method.  Cousteix's  entrainment  relationship  (24)  and 
Lyrio/Fer ziger 's  skin  friction  relationship  (25)  are  used  for 
closure.  The  details  of  the  methods  are  described  in  Ref.  7. 
Tests  of  the  procedure  against  experimental  data  and  against 
other  methods  show  good  agreement  (7). 

The  calculation  provides  the  boundary  layer  displacement 
effect  along  each  of  a  family  of  instantaneous  "streamlines". 
This  is  then  redistributed  to  the  surface  panels  in  the  attached 
flow  regions.  Computed  information  on  the  location  of  separation 
is  assembled  and  passed  over  to  the  wake  routines. 

3.6  BaKfi..  Routines 

When  a  separation  line  is  first  passed  to  the  wake  routines 
a  set  of  wake  panels  is  formed  simply  shedding  the  jump  in  total 
potential  (m  +  $„)  between  each  pair  of  surface  panels  to  either 
side  of  the  separation  line.  At  subsequent  time  steps  all  the 
previous  wake  panels  are  convected  downstream,  each  corner  point 
traveling  for  time,  fit,  along  the  local  computed  velocity  vector. 
Figure  5.  In  this  process  a  new  wake  panel  is  formed  corre¬ 
sponding  to  each  surface  panel  at  the  separation  line  as  in  the 
first  step. 

Several  "housekeeping"  routines  are  necessary  to  complete 
the  wake  model.  First,  the  movement  of  the  separation  line  is 
monitored  from  one  step  to  the  next  and  restricted  if  necessary; 
e.g.,  downstream  movement  is  limited  to  the  distance,  fit  x  local 
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Figure  5.  Wake  Movement  Viewed  in  the  Aircraft  Frame 


transport  velocity.  Also,,  the  movement  of  each  wake  point  is 
examined  for  possible  impingement  into  the  airfoil  surface.  This 
is  an  important  routine.  If  an  impingement  is  detected  and  if 
the  line  if  impingement  is  at  a  low  inclination  to  the  panel 
surface,  then  the  code  simply  "nudges"  the  point  away  from  the 
surface.  If  the  line  of  impingement  is  steep  into  the  panel, 
however,  then  the  wake  is  geometrically  cut  off  at  the  point  of 
intersection,  usually  forming  a  closed  bubble  region.  Doublet 
values  continue  to  be  convected  along  this  boundary  in  subsequent 
time  steps,  and  the  impingement  calculation  is  repeated.  Thus  it 
is  possible  for  the  bubble  to  open  again  under  different  condi¬ 
tions  at  a  later  time;  the  wake  would  then  proceed  to  grow  again. 

When  a  new  wake  configuration  has  been  obtained  the  new  wake 
panels  and  influence  coefficients  are  formed  and  a  new  potential 
flow  solution  is  computed.  Figure  3. 


4 . 0  RESULTS 


As  the  applications  increase  in  complexity,  available  test 
cases  for  comparison  become  scarcer.  It  is  important,  therefore, 
to  have  a  set  of  baseline  test  cases  of  well  proven  results  that 
can  be  checked  from  time  to  time  as  the  method  develops.  Two 
such  cases  are  the  Wagner  function  for  an  impulsively  started 
flat  plate  and  the  Theodorsen  solution  for  an  oscillating  plate. 
These  cases  are  briefly  considered  below  to  establish  the  ac¬ 
curacy  of  baseline  solutions. 

Figure  6  shows  the  growth  of  indicial  lift  and  circulation 
for  a  NACA  0012  impulsively  started  from  rest  at  an  angle  of 
attack  of  0.1  rad.  The  curves  are  compared  with  Wagner's 
function  for  indicial  lift  and  R.T.  Jones'  indicial  circulation 
for  a  flat  plate.  The  calculations,  which  used  31  panels  around 
the  section,  are  in  good  agreement  and  indicate  a  slightly  higher 
trend  which  is  consistent  with  a  higher  steady  state  circulation 
for  the  thickness  case. 

Refinements  in  the  wake  shedding  model  developed  earlier  in 
the  two-dimensional  pilot  code  significantly  reduced  the  com¬ 
puting  requirement  of  these  time-stepping  calculations.  For 
example.  Figure  7  shows  the  effect  on  indicial  lift  of  varying 
the  time-step  size  in  the  Wagner  problem  and  demonstrates  a  rapid 
convergence. 

The  procedure  has  been  tested  also  for  a  NACA  0012  oscil¬ 
lating  in  pitch  about  the  quarter  chord.  The  calculations  com¬ 
pare  favorably  with  the  Theodorsen  flat  plate  function  over  a 
range  of  reduced  frequency.  Figure  8.  The  calculations  were  per¬ 
formed  with  only  16  time  steps  per  cycle.  Figure  9  shows  the 
computed  results,  CL  versus  time,  using  only  4  time  steps  per 
cycle.  This  is  in  remarkably  good  agreement  with  the  16  and  also 
32  time-step/cycle  solutions,  demonstrating  an  extremely  good 
convergence  characteristic. 

Time-stepping  calculations  have  also  been  performed  for 
cases  with  prescribed  extensive  separations.  The  purpose  of 
these  calculations  was  to  check  the  basic  unsteady  circulation 
shedding  model  in  the  potential  flow  code.  For  the  first  set  of 
tests,  the  wake  panel  were  simply  transported  at  the  onset  flow 
velocity  after  the  initial  growth  as  determined  from  the  surface 
conditions  at  separation.  Several  triangular  shapes  were  con¬ 
sidered,  each  starting  impulsively  from  rest  and  proceeding  for¬ 
ward  over  10  time  steps  for  a  total  time  of  t  *  t  U^/h  *  3.0, 
where  h  is  the  triangle  base  height.  Separation  was  prescribed 
at  the  corners.  Figure  10(a)  shows  the  computed  history  of  the 
drag  coefficient  from  pressure  integration  for  a  60  deg  triangle 
with  blunt  face  forward.  A  total  of  40  panels  was  used  to  repre¬ 
sent  the  triangle  surface.  The  calculation  was  repeated  in  the 
presence  of  wind  tunnel  walls  (also  paneled)  with  a  10%  blockage 
ratio.  The  indicated  blockage  correction  is  somewhat  lower  than 
that  given  by  standard  techniques.  Figure  10(b)  compares  the 


IN  HALF  CHORDS 


u 

o 


TRAVELED  IN  HALF  CHORDS 


ISINkt 


lculated  Drag  History  for  a  60°  Wedge  Started 


Calculated  Pressure  Distribution  for  a  60°  Wedge  Started 


computed  pressure  distributions  for  this  triangle  in  and  out  of 
the  tunnel.  The  "base"  pressure  has  only  a  small  variation  and 
is  quite  close  to  experimental  measurements.  Figure  11  shows  a 
summary  of  computed  drag  coefficient  versus  triangle  semi-apex 
angle.  The  calculated  values  are  slightly  high  in  relation  to 
the  experimental  data  collected  from  several  sources  by  Hoerner 
(26)  . 

One  further  case  was  run  for  the  60  deg,  apex-forward  tri¬ 
angle  in  free  air  with  the  full  wake  velocity  calculation  routine 
turned  on  but  without  the  amalgamation  and  redistribution  schemes 
at  this  stage.  The  calculated  CD  for  this  case  falls  below  the 
experimental  value.  Figure  11. 

Although  the  two-dimensional  pilot  program  was  generated 
primarily  as  a  tool  to  examine  the  behavior  of  various  parts  of 
the  dynamic  separation  calculation,  it  has  shown  considerable 
promise  as  a  general  purpose  code  for  two-dimensional  calcula¬ 
tions  involving  separated  flows.  For  example,  it  has  been  ap¬ 
plied  to  the  cross-flow  problem  on  helicopter  booms  and  also  to 
the  download  problems  on  tilt-rotor  aircraft  (27).  Earlier  ex¬ 
amples  (e.g..  Figures  10  and  11)  demonstrated  a  capability  to 
compute  base  pressures  and  drag  coefficients  of  blunt  sections 
using  an  impulsive  start.  An  extension  of  this  to  compute 
spoiler  characteristics  has  also  been  briefly  examined.  Figure 
12  shows  computed  wake  configurations  at  three  steps.  This  is 
for  the  case  of  a  spoiler  deflected  30  deg  on  an  airfoil  at  a  = 
8  deg.  The  final  base  pressure  and  integrated  lift.  Figure  12 
(b)  and  (c),  respectively,  are  in  good  agreement  with  experimen¬ 
tal  measurements  (28).  The  calculated  values  represent  an 
average  value  over  the  last  few  time  steps  as  the  solution  had 
started  to  oscillate.  The  amplitude  on  Cl  is  about  0.1,  but  the 
calculation  ought  to  be  continued  for  a  longer  time  to  examine 
whether  a  pattern  between  upper  and  lower  vortex  formation  is 
established.  This  application  could  be  extended  further  to  ex¬ 
amine  pitch  rates  and,  with  a  straightforward  extension  of  the 
code,  rates  of  spoiler  deployment.  Such  an  extension,  involving 
relative  motion  between  parts  of  the  configuration,  would  also 
allow  treatment  of  pitching  airfoils  between  channel  walls  to 
assess  the  effects  of  unsteady  blockage  in  pitch-up  or  oscilla¬ 
tory  experiments. 

An  experimental  data  case  for  an  airfoil  oscillating  in 
pitch  (29)  was  run  and  the  computed  lift  variation  with  a  com¬ 
pared  with  the  measured  data  in  Figure  13.  The  airfoil  is  a  NACA 
0012  and  is  oscillating  in  pitch  about  the  quarter-chord  line 
with  a  =  8.1  deg  +  4.9  deg  sin  (0.2t);  4*6.,  below  the  dynamic 
stall  onset.  Reynolds  number  was  4  x  10  .  This  reduced  frequen¬ 
cy  condition  is  very  close  to  the  changeover  from  a  lead  to  a  lag 
situation  and  so  there  is  only  a  small  difference  between  the 
upswing  and  downswing  curves. 
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REDUCED  FREQUENCY:  0.2 

Re:  4  x  106 
M„:  0.3 

a:  8.1*  ♦  4.9*  sin  «<ft 
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Figure  13.  Comparison  of  Calculated  and  Measured  Lift  on 
a  NACA  0012  Airfoil  Oscillating  in  Pitch  about 
the  Quarter  Chord. 


A  test  calculation  was  performed  for  a  NACA  0012  section  in 
a  state  of  pitch-up  from  10  deg  to  30  deg  at  a  rate#  ac/2Uoo  * 
0.175.  The  surface  was  represented  by  30  panels.  Separation 
point  movement  on  the  upper  surface  was  prescribed  based  on 
boundary  layer  calculations  performed  on  the  pitching  case  with 
just  the  trailing-edge  wake.  (This  case  was  run  prior  to  the 
automatic  coupling  between  the  viscous  and  inviscid  codes.)  The 
motion  was  started  impulsively  from  rest  and  proceeded  through  10 
time  steps. 

Figure  14(a)  shows  a  sample  of  the  computed  wake  shapes  and 
demonstrates  a  reasonable  numerical  behavior.  Sample  pressure 
distributions  are  shown  in  Figure  14(b).  The  passage  of  the 
leading-edge  vortex  is  clearly  shown.  This  is  associated  with  a 
local  region  of  reversed  flow,  Figure  14(b)  (ii).  These  were  pre¬ 
liminary  test  calculations  aimed  at  exploring  the  numerical  be¬ 
havior  of  the  calculation  procedure  and  potential  flow  model. 

Experimental  measurements  of  airfoils  undergoing  constant 
rate  pitch-up  motions  have  been  taken  at  the  Frank  J.  Seiler 
Research  Laboratory  (19).  Three  cases  are  considered  briefly 
here  for  correlation  purposes.  In  these  cases  the  airfoil  is 
pitched  up  from  a  =  0  to  approximately  1  radian  at  a  constant 
pitch  rate  and  then  held  at  constant  angle  of  attack.  Three 
pitch  rates  are  considered  with  normalized  pitch  rates,  k 
of  0.047,  0.089  and  0.133.  Figures  15(a),  (b)  and  (c), 
respectively,  show  the  comparison  between  calculated  and  measured 
CL  ~  a  and  CD  “  a  characteristics  for  these  cases.  For  the  low 
pitch  rate,  k  »  0.047,  the  comparison  is  very  good  up  to  about  30 
deg,  but  then  deteriorates.  The  calculated  CL  remains  fairly 
constant  with  a  until  the  pitch  rate  drops  to  zero  while  the 
experimental  curve  falls  markedly.  The  calculated  rise  in  drag 
with  a  has  a  steady  rate  in  the  30  deg  to  60  deg  range  while  the 
experimental  measurements  include  a  substantial  increment  above 
this  rate  peaking  at  about  a  *=  40  deg.  There  is  a  good  agreement 
between  calculated  and  measured  lift  and  drag  value  for  the  final 
"steady  state"  conditions  at  a  =  60  deg. 

One  possible  reason  for  the  departure  of  the  calculated  lift 
and  drag  in  the  latter  part  of  the  pitch-up  phase  is  the  modeling 
of  the  leading-edge  vortex  roll-up.  The  amalgamation  and  redis¬ 
tribution  schemes  that  were  installed  to  stabilize  the  dynamic 
wake  calculations  were  not  performing  in  a  consistent  manner. 
The  vortex  core  on  the  upper  sheet,  therefore,  did  not  "condense" 
early  enough  in  the  amalgamation  procedure;  consequently,  the 
vortex  formed  just  downstream  of  the  trailing  edge  and  did  not 
closely  interact  with  the  airfoil  surface  during  the  pitch-up 
phase.  This  tendency  was  still  present  at  the  higher  pitch-up 
rates,  k  *  0.089  and  k  =  0.13,  but  to  a  lesser  extent.  In  these 
cases  the  calculated  lift  and  drag  characteristics.  Figures  15(b) 
and  (c),  respectively,  are  in  very  good  agreement  with  experi¬ 
ment.  The  tendency  for  the  measured  lift  to  peak  at  about  a  *  30 
deg  is  also  shown  in  the  calculated  results.  These  calculations 
were  not  continued  at  amax  for  a  sufficient  time  to  enable 
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Figure  14(a).  NACA  0012  Starting  Impulsively  from  Rest  Pitching  from  10°  to  30°: 
<Sc/2U  =  0.175. 
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Figure  14(b)(i).  NACA  0012  Starting  Impulsively  from  Rest  Pitching 
from  10°  to  30c :  uc/2U  =  0.175;  t  =  0.6. 
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Figure  15.  Comparison  of  Calculated  and  Measured  Lift  and 
Pressure  Drag  on  a  NACA  0012  Section  During 
Pitch-up  Motion  about  x/c  -  .317. 


"steady  state"  conditions  to  be  reached.  As  stated  earlier, 
because  of  problems  with  the  amalgamation  routine,  a  strong 
vortex  core  did  not  "condense"  for  these  cases  and  so  the  upper 
surface  suction  peak  seen  in  the  experimental  measurements  (19) 
did  not  materialize;  rather,  a  smeared  suction  peak  appeared 
because  of  the  more  diffuse  region  of  shed  vorticity.  Conse¬ 
quently,  although  the  integrated  lift  and  drag  are  in  good  agree¬ 
ment,  the  pitching  moment  characteristic  (not  shown)  is  not 
satisfactory  at  this  time.  This  upper  surface  suction  peak, 
which  is  associated  with  a  reversed  flow  region  under  the  vortex, 
was,  in  fact,  computed  in  earlier  preliminary  calculations  (see 
Figure  14)  involving  a  high  pitch  rate,  k  =  0.175.  In  this  case 
a  vortex  core  condensed  early  in  the  calculation  (Figure  14(a)). 

Earlier  tests  of  the  pitch-up  case  were  repeated  with  dif¬ 
ferent  conditions  to  compare  with  more  recent  data  by  Helin  and 
Walker  (20)  at  the  Seiler  Laboratories.  Figure  16(a)  shows  a 
comparison  of  wake  geometry  from  two  a  conditions  straddling  the 
a  =  45  deg  condition  of  a  flow  visualization  photograph.  The  wake 
"boundary"  in  the  calculation  shows  good  qualitative  agreement 
with  the  photograph.  The  pitch  rate,  itc/lV^  **  0.2,  and  the  pitch 
axis  is  the  quarter-chord  line.  For  the  same  case  a  history  of 
the  velocity  magnitude  at  an  upper-surface  point  at  x/c  «  0.7  is 
shown  in  Figure  16(b)  as  a  function  of  a.  This  is  in  remarkably 
good  agreement  with  the  measurement. 

While  the  two-dimensional  program  has  been  used  to  examine 
and  develop  the  various  routines  required  for  the  coupled  dynamic 
separation  calculations,  the  three-dimensional  pilot  code  de¬ 
velopment  has  been  following  closely  behind.  The  unsteady 
boundary  layer  calculation — which  is  performed  along  computed 
surface  (external)  streamlines  at  each  time  step — has  been  fully 
coupled  with  the  unsteady  inviscid  program.  Test  cases  have  been 
performed  and  compared  with  experimental  data  from  the  DFVLR-AVA 
in  Gottingen.  These  experiments  were  conducted  as  part  of  a 
cooperative  agreement  between  the  DFVLR  Institute  of  Aero- 
elasticity/West  Germany  and  NASA  Langley  Research  Center.  Figure 
17  compares  the  calculated  and  measured  real  and  imaginary  pres¬ 
sure  distributions  at  a  70%  spanwise  station  on  an  AR  *  4  rec¬ 
tangular  wing  undergoing  pitch  oscillation  about  the  quarter 
chord  a  *  7.9  deg  +  1.0  deg  sin  (.2t).  (Reynolds  number  is 
1.35  x  10  .)  The  potential  flow  solution  is  also  shown  to  indi¬ 
cate  the  extent  of  the  viscous  correction.  The  complete  solution 
is  in  very  good  agreement  with  the  measurements.  This  is  still 
true  for  the  condition,  a  =  12  deg  +  1.0  deg  sin  (.3t),  which  is 
approaching  the  condition  of  dynamic  stall  onset;  a  pressure 
deviation  is  apparent  near  the  leading  edge. 


(a)  Wake  Geometry  at  Two  Time  Steps 


(b)  Variation  of  Velocity  at  x/c  =  0.7  Versus  Angle  of  Attack 
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Figure  17(a).  Comparison  of  Chordwise  Pressure  Oistribut 
between  Computed  and  DFVLR  Data  (uQ  =  7.9* 
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5.0  CONCLUSIONS 


A  system  of  routines  has  been  developed  coupling  an  unsteady 
time-stepping  potential  flow  panel  code  with  unsteady  integral 
boundary  layer  methods.  The  routines  include  treatment  of  the 
growth  of  multiple  vortex  sheets  representing  free  shear  layers 
in  the  wake.  The  location  of  separation/  predicted  by  the 
boundary  layer  code,  is  allowed  to  move  with  time.  The  procedure 
is  completely  assembled  in  a  two-dimensional  pilot  code  and  shows 
encouraging  results  for  pitch-up  conditions  beyond  the  dynamic 
stall.  Unsteady  viscous/inviscid  calculations  in  the  three- 
dimensional  pilot  code  are  in  close  agreement  with  experimental 
measurements  up  to  dynamic  stall  onset. 
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